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Abstract 

By means of an effective potential associated to a constrained equilibrium measure and apt to study 
frozen systems, we investigate glassy freezing in simple liquids in the hypernetted chain (HNC) approxi- 
mation. Differently from other classical approximations of liquid theory, freezing is naturally embedded 
in the HNC approximation. We get a detailed description of the freezing transition that is analogous to 
the one got in a large class of mean-field long range spin glass. We compare our findings with Monte 
Carlo simulations of the same system and conclude that many of the qualitative features of the transition 
are captured by the approximated theory. 



1 Introduction 

Cooled at fixed rate, supercooled liquids stop flowing on observable time scales when the glassy transition 
temperature T g is met JTJ. At that cooling rate dependent temperature these systems get out of equilibrium. 
The large scale motion of the molecules is frozen, and consequently the entropy (in the sense of the logarithm 
of the phase space accessible to the system) is discontinuously reduced. This discontinuity is equal to 
the configurational entropy at temperature T g and is supposed to vanish if the system could be kept in 
equilibrium down to the point of the ideal glass transition To. This would be observable only for infinitely 
slow cooling, and would be a real thermodynamic transition. 

At temperatures smaller then T g these systems finds themselves in a region of the configuration space 
having vanishing weight in the Boltzmann distribution^, the values of the extensive quantities being far from 
these at equilibrium. One can also expect that they will remain conformationally close to the configuration 
' y reached at T g where flow stopped. On the other hand small scale vibrations of the atoms are free to 
thermalize at the actual temperature T of the thermal bath. This situation, with extreme separation of 
time scales would be naturally described by a statistical ensemble where the slow degrees of freedom are 
quenched, and the fast degrees of freedom thermalize at temperature T. One can then define a conditional 
statistical ensemble, as a Boltzmann-Gibbs measure for fixed distance from the point y. The free-energy 
associated to this distribution is a function of the constrained distance, and is a natural "effective potential" 
for glassy systems. This effective potential can be computed with the replica method even in non disordered 
models, the role of the quenched variables being played by the reference configuration y. This reduces the 

1 Strictly speaking the in the metastable region the liquid has zero weight, the Boltzmann distribution being concentrated 
on crystal configurations. In this paper we simply neglect the existence of the crystal and we imagine that in the supercooled 
region the measure is concentrated on liquid configurations. This situation can also be realized by not considering in the 
partition sum the crystal like configurations or by modifying the potential in such a way that the crystal get an high free 
energy. 



problem to the study of the free-energy of a multicomponent mixture, in which an analytic continuation on 
the number of components has to be performed at the end of the computation. 

In previous work j|, [3|, |j, it was studied the shape of the potential, and the implications for the glass 
transition in long-range spin glasses, and supported the generality of the picture in numerical simulations 
of a binary mixture model Q . In this paper we extend the analysis to models of simple liquids in the HNC 
approximation, and show how this approximation, devised to study the liquid phase naturally describe 
glassy freezing || . The same would not be true for other classical approximations of liquid theory as the 
Percus-Yevick or the Mean Spherical approximations. The implementation of the replica formalism for the 
effective potential in the HNC approximation is similar to the one used by Given and Stell §1 to study 
liquids in random quenched matrices. It also bear resemblance with the one used by Zippelius and coworkers 
to implement random crosslinking in models of vulcanization. The main difference is that, while in these 
cases the replica method is used to deal with external quenched disorder, in our case we use quenched 
degrees of freedom to probe the configuration space of systems that freeze even in absence of quenched 
disorder. 

The approach of this paper is complementary to the one put forward recently in ||. There it was 
shown how, combining HNC and replicas, one could reveal the glassy transition and find the properties 
of the system below the glass temperature. In this paper we will discuss the nature of the freezing in the 
HNC approximation finding a scenario very similar to that of mean-field spin glasses. The approximation 
is constructed in such a way that the critical density automatically coincides with the one obtained in (5j. 
The advantage of the "effective potential" framework with respect to that of ^ is to make conceptually 
clear the introduction of the replicas in the theory and to make testable predictions on the behaviour of 
the system at density less the critical one (as we shall see below) when we introduce a potential among two 
copies of the system. On the other hand we will see that for high values of the density the simpler approach 
of this paper, where we neglect replica symmetry breaking, leads to inconsistency, and there one need to 
resort to the approach of j|] for a coherent theory. 

The explicit computations are performed for the hard sphere potential. We support our findings with 
Monte Carlo simulations of the same system. A short account of our result has appeared in 0. 

We organize the paper as follows: in section 2 we discuss the construction of the effective potential and 
we briefly review the results obtained for long-range spin glasses. In section 3 we discuss the potential for 
simple liquids in the HNC approximation. Section 4 is devoted to the presentation of the theoretical results 
on the hard sphere system, that in section 5 we compare with the numerical simulations. Finally in section 
6 we present some conclusions and perspectives. 

2 The effective potential 

In this section we review the construction of the effective potential @, ||[ f|]. For definiteness we discuss 
the case of a simple liquid composed by N identical point-like particles in a volume V, described by their 
coordinates x = (xi, xn), and interacting via a pair potential 4>(xi — xj). Suppose that, undergoing a 
cooling process from the liquid phase, the system falls out of equilibrium at a temperature T g and remains 
stuck in a region of the configuration space having vanishing weight in the Boltzmann-Gibbs measure. This 
commonly happens at the glassy transition of supercooled liquids, where the liquid stops flowing: large 
scale motion is frozen, while small scale motion of the atoms (vibration) can still equilibrate even below T g . 
In these conditions the observed values of extensive quantities can be far from their canonical equilibrium 
values, while keeping the external parameter constant, they do not vary over the laboratory time scale. 
It is appropriate then to restrict the measure in configuration space to the vicinity of the configuration y 
reached when crossing T g . 

In order to do that we need to define a notion of similarity (or codistance) among configurations, q(x, y), 
that, with reference to spin glass terminology, we call overlap. The appropriate definition of the overlap 
depends on the problem at hand, and it has to be such that to similar configurations correspond high values 
of q (with normalization q(x,x) = 1) and to very different configurations values close to zero. 

In our particle system an appropriate definition can be 

q{x,y) = -fiY, w (\ x i -Vj\) (!) 

i,3 



where w(r) is a function close to one for r < err® and close to zero for r > aro, with tq being the radius 
of the particles and a a number e.g. of the order of 0.3, such that couples of particles at small distances 
in the two configurations contribute positively to q. Specifically in the unit radius hard sphere problem of 
section 4 we will use w(r) = 9{r — 0.3). As q(x,x) = 1 and q(x,y) < 1 we can define a sort of distance as 
d(x,y) = 1 — q(x,y). In the following we will speak indifferently about the two quantities, remembering 
that high overlap means small distance and vice-versa. 

Having now the definition of q(x, y) we can define a restricted Boltzmann-Gibbs distribution as 

P{x\y) = exp (-0H(x)) 5(q(x, y) - q) (2) 

where Z((3,y) is the integral over x of the numerator of (|2|). 
Three comments are in order. 

• The value of q that appears in (||) is at this stage arbitrary. However, the system at temperature 
T will tend to adjust itself and select a given natural distance from the configuration y, according 
to local free-energy minimization. The selection of q can be well understood in a mean-field picture, 
and has been discussed in the case of long range spin glasses in . At low temperature one expects 
metastable states in configuration space. This corresponds to a two-minima structure of the effective 
potential, with one minimum at low q representing the typical overlap among configurations belonging 
to different metastable states, and one minimum at high q representing the typical overlap among 
configurations in the same metastable state. This last is the q that would be naturally chosen by the 
system. 

• The second comment concerns the dependence of the measure (|2|) on the reference configuration y. At 
a first sight, as different cooling experiment would produce different configuration y, it would appear 
that the measure (^) could be hardly of any use. However, y is supposed to be a configuration typical 
with respect to the Boltzmann-Gibbs probability at temperature T g , n(y) = exp(—f3 g H(y))/Z((3 g ), 
and we can expect the extensive quantities computed from (|2|) to be self-averaging (i.e. y independent) 
in the thermodynamic limit. The role of the configuration y is analogous in this construction to the one 
of the quenched variables in disordered systems. In this sense, the measure (0) is an implementation 



of the idea of "self-generated disorder" often advocated for structural glasses P, 12]. 



• The third comment concerns the selection of the temperature T g . For that we do not have any a-priori 
criterion, as it is a quantity that in experiments depends on the cooling rate. We have thought to 
the temperature of the configuration y, that we will call T' in the following, as the glass transition 
temperature in the purpose of illustrating the physical situation that we have in mind. As the matter 
of fact our construction is well defined for arbitrary T', and interesting results are obtained even 
for T' = T. In this paper we will limit out analysis to this case using the measure (|2]) as powerful 
probe of configuration space. We insist however on the conceptual importance of considering two 
temperatures, and we will often refer to results obtained in spin glasses for this more complicated 
case. In the previous discussion and part of the following we have considered the temperature as the 
only external parameter. It is clear that mutatis mutandis analogous considerations hold for any other 
control parameter, as the density of section 5. 

The object on which we will concentrate our attention is the free-energy associated to the distribution 
@ that, invoking the self-averaging property we can write as: 

V(q,P,P') = -§^7) / dyexp(-f3'H(y))lo g y dx exp(-f3H(x))5(q(x, y) - </)} (3) 

As the constraint implied by the delta function is global, we can enforce it through a Lagrange multiplier; 
and considering the quantity 

F(e,P,P') = -^^ry / dy eM-P'H(y)) log {J dxexp(-/3[H(x) - eq(x,y)])} (4) 



we find that V and F are related by the Legendre transform 



V(q,P,P') = mm(F(e,(3,P') + eq). (5) 

For practical purposes it is more convenient to work with F than with V, while the data are more easily 
interpreted in terms of V . We will pass freely from one representation to the other in the following. 

In order to deal with the average of the logarithm in (Q) we resort to the replica method. This consists in 
evaluating the moments Z r for integer r and computing average the logarithm from an analytic continuation 
to non integer r from the formula log Z = lim r _>o Z ' ? r 1 ■ Explicitly we can write: 

W=-T f dx Qe M-P'H{x Q ))/Z{(i') { dx^.Axrexp (-f3[J2H(x a ) - eJ2Q^a,x ))) . (6) 

J J V o=l a=l ) 

where we have written y = xq. The problem is reduced to the computation of thermodynamics for a mixture 
of r + 1 components in the limit r — > 0. Notice the non-symmetric role played by the replica xq and the 
replicas x a for a > 1. This implies that while there is symmetry under permutation of replicas with positive 
index there is not symmetry under interchange of the replica xq with the others. Although we have r + 1 
replicas, the symmetry of the problem is only S r , becoming S r +i only for e — > 0. Technically, our approach 
is similar to the one of Goldbart and Zippelius [13| to study vulcanization of rubber and the one of Stell et 



al. for liquids in random quenched matrices, where also one finds a number of replicas that tends to one. 
However, in their case real quenched disorder is present, while in our case we use the auxiliary configuration 
y to restrict the Boltzmann-Gibbs measure to small regions of configuration space. 

Before discussing in the next section the application of the present formalism to simple liquids in the 
hypernetted chain approximation, let us describe briefly what one can expect for the effective potential 
when a glassy transition occurs, considering for simplicity the case T' = T. In the supercooled phase the 
diffusion constant becomes lower and lower as the temperature is lowered. The "cage effect" takes place: 
the molecules get trapped for long times before they can diffuse. When the glass transition is met diffusion 
is completely stopped, at least on human time scale. The entropy associated to diffusion is lost at the 
transition. Ergodicity is broken and the configuration space is effectively split into an exponentially large 
number of practically mutually inaccessible regions. 

Making the approximation that the time to jump out of these regions is infinite, it is natural to expect 
the effective potential (|2|) to have two minima. One corresponding to the typical (low) overlap among 
configurations belonging to different regions in configuration space, and another corresponding to the typical 
overlap (high) of different configurations belonging to the same region. The number of these regions, or 
metastable states J\f is related to the configurational entropy £ by the relation: J\f = exp(iV£). The 
probability of x to be in the same metastable state of y will be in such conditions 1/Af = exp(— NT,). 
Consequently the relative height of the high q minimum with respect to the low q one has to be equal to 
T£. This picture is realized and it has been discussed in ||| in a large class of long range spin glass 
model. The analysis of these models tells us that the picture has to be refined a little. To each disconnected 
region one can associate a free-energy /, with an energetic part and an entropic part. Defining £(T, /) the 
logarithm of the number of these regions as a function of /, one finds that, at low enough temperature, £ 
is different from zero in a finite temperature dependent interval I(T) = [f m (T), fu (T)]. The states that 
dominate the partition function at temperature T are such if the quantity 

F = f — T£(T, /) (7) 

is minimum jnj. The study of the effective potential for T ^ T' shows that individual states do not 
disappear when the temperature is changed, but they remain stable for large ranges of temperatures f^, fig , 
|l6[| . As the temperature is lowered, states with lower and lower £ are selected in (0) until, for a temperature 
T s with £ = are reached and the partition starts to be dominated by the lowest states. 



3 The HNC approach 

Let us now start the discussion of the implementation of the HNC approximation in the approach outlined 
in the previous section. As we stressed the use of the replica method reduces the problem of the evaluation 



A 



of the effective potential to the one of an r + 1 component liquid mixture, in which there is a privileged 
component with which all the other replicas interact via the potential 



n n 1,N 

-NeY / q(x a ,x ) = -eY / Y, w ( x °- x j) ( 8 ) 

a=l a=l i.j 

The basic quantities of the theory are the pair correlation functions among replicas 

PaPb9ab{x, y) + p a 8ob8{x -y) = ^5(x? - x)5(x) - y) (9) 

A partial resummation of the Mayer expansion allows to write a self-consistent expression for the free-energy 



- 2(3F H nc = / d d x ^2 p a Pb9ab{x) [\ogg a b{x) - 1 + P a <t>(x)5 ab ] 

a,b=0 
r 

+2(5e PQpagva(x)w{x) + Tr L(ph) (10) 
a=l 

where h ab = g ab — 1 is the connected correlation function and L is an operator in physical and replica space, 
defined by 

L(u) =u-u 2 /2-\og{l + u). (11) 

We have also put (3q = /?', po = p' and f3 a = (3, p a = p for a > 1. The free-energy (|T~0| ) has to be extremized 
with respect to the Safe's and the terms of order r have to be extracted. The extremum conditions can be 
cast in the form: 

9ab{x, y) = exp(-/3 a (j)(x, y)5 ab + {5 0a (l - 8 ab )(5 b + 8 0b (l - 5 0a )(3 a } ew(x) + h ab (x) - c ab (x)) (12) 
with g and c related by the Ornstein-Zernike relation 

hab(x) = c ab (x) + ^2 dy h ac {x - y)p c Ccb{y) (13) 
c =cr 

As usual in the replica method one needs a parameterization of the matrix g ab that allows the analytic 
continuation to r — > 0. On the basis of symmetry considerations analogous to the one of 0, one can propose 
the structure: 

{5oo a = b = 

gw a = 0, b + or b = 0, a + (14) 

The r x r sub-matrix g* b can be either replica symmetric or have an ultrametric structure |f[G| ], in the 
following we will limit ourselves to the replica symmetric structure 

j 9n a = b , . 

9ab= \9i2 a + b. ^ 

which coincides with the choice done by Given and Stell ||. In the case of the p-spin model, the replica 
symmetric choice gives the correct result for the effective potential in the high and in the low q regions and 
in particular around the minima. However it was found an intermediate q region where "one step replica 
symmetry breaking" in the matrix g* b was necessary to compute correctly the effective potential. Although 
we do not explore here the possibility of solutions with a structure more complicated then ( |l5| ) and we limit 
ourself to the study of the high and low q parts of the effective potential, we warn the reader that replica 
symmetry breaking has also to be expected in this case for intermediate q. 

The interpretation of the different elements of the g ab matrix with the replica symmetric ansatz is 
straightforward. The element 500 represents the pair correlation function of the free system; as such the 
equation determining it decouples from the other components in the limit r — > 0. In turn, gn represents the 



pair correlation function of the coupled system. gio is the pair correlation among the quenched configuration 
and the annealed one, while g\ 2 represents the correlation between two systems coupled with the same 
quenched system. This is the analogous of the Edwards-Anderson order parameter in disordered systems, 
and represents the long time limit of the time dependent autocorrelation function at equilibrium Q. 

Straightforward algebra shows that the equations (|i~3|) reduce to the ones proposed in || . In particular 
one finds that, as it should, the equation for goo, describing the correlation function of the quenched replica 
decouples from the other and coincides with the usual HNC equation in absence of replicas. 
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+ Pi{h u (x - y)c 12 (y) + h 12 (x - y)cu(y) - h 12 (x - y)c 12 (y)(2 - r))) 
The overlap can be expressed in terms of the correlation function yio, and reads 

Q = PJ dxw(x)gio(x) = Anp r 2 w(r)g w (r). (17) 

According to the discussion of the previous section, we will associate glassy behavior to non convexity of 
the function V(q), and in particular to the existence of multiple solutions q(e) for e — > 0. In the liquid phase 
we can expect instead q{e) to be a single value function and a convex effective potential with a minimum 
at q = qo = p J dxw(x), corresponding to the absence of any structure in yio, i.e. gio(x) = 1 for all x. A 
strong coupling e, attracting the system towards the configuration y will force a structure in the yio, which 
will have a higher peak in x = the higher is the coupling. Similarly g± 2 will acquire a structure: if two 
different systems are similar to the configuration y they will also be similar to each other. When the system 
freezes, it will exist a solution to the HNC equations in which gio and gi 2 will have a structure even for 
small and vanishing e. 

The value of e at which we find a solution with non-zero gi 2 coincides with dynamical critical density 
of and it correspond to the phase transition point in a mode-coupling approach [|| . Indeed in this mean 
field approach ergodicity starts to be broken exactly at this point. 



(16) 



4 Results for HNC hard spheres 

In this section we discuss the picture coming from the integration of the equations ( |l2|p^ ) in three dimension. 
We present systematic data in the case of the hard sphere potential 

0( r ) = I n r < I an d with w(r) = < n T < " (18) 
^ w 1 r > 1 w 1 r > 0.3 v ; 

The hard sphere potential has been chosen for practical convenience, the glassy transition picture that 
will emerge can be strongly expected to be very general. We have verified in non systematic investigations 
that the same picture indeed holds for soft sphere systems with (j>(r) = r~ 12 . The value 0.3 that appears in 
the definition of w has obviously nothing fundamental, and we have checked that the picture is insensitive 
to its precise value. The hard sphere model has no temperature, and the control parameter is the density. 
Numerical work report a glassy phase for values of the density higher then 1.15. 

We have solved the saddle point equation (12,1^) by iteration for various values of the density and the 



coupling. For fixed density we start the integration of the equation at low (respectively high) coupling e 
where we know the solution and we increase (respectively decrease) it at small steps. In this way we can find 
the curves of q as a function of e, and reconstruct from (|io|) the effective potential V(q). We checked that 
up to a constant V(q) = J q dq'e(q'). In the low density region we were able in this way to fully reconstruct 
the shape of the potential. For higher densities, we could just reconstruct in this way the high and the low 
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Figure 1: The behavior of q as a function of e for HNC hard spheres for p = 1.14, 1.17, 1.19, 1.20. For high 
enough density q is a multivalued function of e. We have shown only a portion of the curve in the region 
where it is multivalued. For graphical transparency in this and the next figure we have joined with a line 
the branches corresponding to the same density. 



parts of the effective potential. This is however enough to get a fully detailed picture of the freezing in the 
system. In figure [l] we present the curves of q as a function of e for various values of the density, while in 
figure [2] we plot the corresponding curves V{q). At low density e is a monotonic function of q, testifying 
ergodic behavior of the system. The potential V(q) is convex and has a single minimum for e = 0, where 
the value of the overlap is qo = /947r(0.3) 3 /3 = p x 0.113, corresponding to 510(2;) = 1 for all x. 

Interesting behavior appears for densities higher or equal to p cr ~ 1.14. At p cr the function q(e) begins 
to be multivalued, the potential looses the convexity property and a phase transition among a low q and 
a high q phase can be induced by a coupling. The point (/9 cr ,e cr ), with e cr = 0.305 is a critical point of 
second order phase transition, from which it departs a first order phase transition line ej r (T) (fig. The 
term —Neq(x, y) in the Hamiltonian implies an energetic advantage for the configurations x close to y and 
induces a transition between a high q "confined" phase with high energy and a low q "deconfined" phase 
with high entropy. The transition line tells that generic equilibrium configurations lie in metastable states 
for densities higher that p c . For p cr < p < p c = 1.17 the metastable states have a finite life, and a coupling 
e > ctr(T) is needed to stabilize them. At p c a minimum develops in the potential, and the metastable 
states have an infinite time life. The equation q = q(e = 0) has more then one solution. It has been shown 
by explicit calculation in Q that a second minimum in the effective potential implies that in the equilibrium 
dynamics the system remains confined in a region with a large overlap with the initial state. In figure § we 
plot the function gio(r) for p = 1.20 and various values of e in the high q and the low q solutions. We see 
how the low q solution has little structure (gio ~ 1) while the high q solution has a very pronounced peaks 
for integer values of r. 

In ordinary cases multiple minima in the effective potential as a function the order parameter signal the 
presence of different (stable or metastable) phases with different qualitative characteristics (e.g. liquid and 
gas). Here the implications of the two minima structure are different. The appearance of the secondary 
minimum signals the breaking of the ergodicity, i.e. the split of the support of the Boltzmann-Gibbs measure 
into many, mutually inaccessible, regions. It is easy to realize the link among the two minima structure and 
such non-ergodic situation. If we suppose that the different region have typically all the same distance, one 
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Figure 2: The effective potential for HNC hard spheres. From to to bottom p = 1.0, 1.14, 1.17, 1.19, 1.20. 
For low density, high up in the liquid phase the potential is convex. In the glass phase two minima are 
present. 
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Figure 3: Phase diagram in the plane e — p. A first order transition line terminating in a critical point 
separates a low q from a high q phase. 
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Figure 4: The function gio(r) for p = 1.20 and e = 0.10, 0.15, 0.20. The three curves on the top correspond 
to the high q solution, the three curves on the bottom to the low q one. 



shall have a minimum corresponding to that distance, and another minimum corresponding to the typical 
distance among configurations in the same region. In this perspective the minima are different manifestation 
of the same phase. Coherently with this picture, the internal energy in the two minima should be the same. 
Obviously this last sentence has no meaning for hard spheres where the internal energy is not defined, but 
it can be easily checked in soft sphere systems. The difference of height between the two minima, AV can 
also be understood in this perspective as being due to the fact that the number of regions in which the 
configuration space has split (7V)is exponentially large in the number of particles M = exp(iVE), and each 
region carry a vanishing weight in the measure. In that conditions, chosen y in a region at random, the 
probability that x falls in the same region, which should be exp {-PAV) is equal to = exp(-iVE). 
We find then that E, to be identified with the configurational entropy, is related to AV byQ 



AV = T£. 



(19) 



We have then a method to compute the configurational entropy, that we plot in figure @. An equivalent 



method has been proposed in 2C ] . We see that £ is a decreasing function of the density and vanishes at a 
density p s ~ 1.203,^] according the scenario of Gibbs and Di Marzio of the glass transition, and analogously 
to long range spin glasses. At each value of the density one choose these states such that the total balance 
between / the internal free-energy of the region and the configurational entropy is such to minimize the 
total free-energy. 

Above p s the equations we are considering give the clearly unacceptable result of a negative configura- 
tional entropy, and the approach must be modified. Previous experience in spin-glasses tells us that the 
paradoxical behavior has to be ascribed to an incorrect description of the quenched replica y. For p > p s 
the lowest / states are chosen. These carry a finite Boltzmann weight and a correct description has to take 
this into account. In the correct HNC approach to the high density regime the quenched replica should 
be described by the replica formulation of Mezard and Parisi with replica symmetry breaking H. The 



2 the previous considerations should be modified for T' 7^ T or p' 7^ p. In that case the secondary minimum reflects the 
properties of the states of equilibrium at the primed values of the parameters when "followed" at the non primed values. 

3 The values of p c and p a are compatible with those found in [S], indeed the potential method reproduces the results of 
replica symmetry breaking for the static and the dynamic critical density. 
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Figure 5: The configurational entropy £ as a function of p. 

HNC approximation, originally devised to study liquids at equilibrium, naturally embeds glassy behavior 
in a glassy transition scenario completely analogous to the one of disordered models with "one step replica 
symmetry breaking" . 

In many senses the picture is genuinely mean-field like. In fact, it has been stressed many times that in 
real systems metastable states with infinite life time do not exist, and a mechanism should restore ergodicity 
between p c and p s . Moreover it is easy to realize that the effective potential we have defined must be convex 
beyond mean field. This can be seen constructing configurations with overlap inhomogeneous in space and 
with a free energy lower or equal to the convex envelope of the potential. We expect however a reflex of 
the mean- field structure in real systems. What it is seen as a sharp transition in mean field can still be 
observed as a crossover in real life and some of the prediction of mean field can be expected to hold in finite 
dimensional systems. In particular the existence of a first order transition line in the plane p — e, which 
depends on the existence of metastable states, regardless if their life is finite or infinite we expect to hold 
in real (or realistic) systems. In next section we will submit that to test. 

Before leaving this section let us remark that, how we show in some detail in the appendix, the possibility 
of a glass transition, associated to nontrivial gw(x) for e — > is excluded by other classical approximations 
of liquid theory, the Percus-Yevick approximation and the Mean Spherical Approximation. 

5 Numerical Simulations 

In order to test the predictions of the transition scenario of the previous section we have performed Monte 
Carlo simulations of a system of hard spheres in three dimension coupled with a quenched configuration. 
We have done simulations with a number of particles N ranging from 256 to 1024 and we have not observed 
any significant dependence on the volume. 

To generate the quenched equilibrium configurations at fixed density we start of N particles of zero 
radius in a box with periodic boundary conditions, and we let the radii grow until two particle do not get in 
contact. At that point we make a Monte Carlo sweep, i.e. we move the particles of random amount and we 
accept the change if two spheres do not overlap; the size of the proposed move is fixed is such a way to have 
0.4 average acceptance. We iterate the procedure until the desired density is reached. The volume and the 
radius (r) are at the point rescaled in order to have r = 1. We thermalize then the system for 4000 Monte 

1 n 




Figure 6: The behavior of q as a function of e for a system of 258 particles and p = 1.04. The different curves 
correspond different thermalization times 2 k for each value of e. ^From top to bottom k = 17, 19,21,23. 
For larger thermalization times the system seem to develop a first order jump in q. 
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Figure 7: The same as figure |6] but with an higher density, p = 1.10 and N = 1024. The values of k are 
k = 17, 19, 21. Notice that for k = 17 the system has not relaxed to the low q value even for e = 0.1. 
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Figure 8: The same as figure ^ with p = 1.14 and N = 256. The values of k are k = 17,19,21. As 
expected, for higher densities the transition is pushed to lower values of e 



Carlo sweeps and use the configuration y reached as external field for our coupled replicas experiment. The 
relatively short thermalization is chosen in order to avoid crystallization. We have careful verified that the 
equilibrium configuration does not have any signs of crystallization by monitoring his correlation function, 
which has a smooth minimum around 1.4, has it should do in the liquid phase. 

Having generated the configuration y we start the evolution of a coupled system x. For various densities, 
we start the evolution from the configuration y with an high value of e and decrease the value of e in units of 
5e, making 2 k Monte Carlo iterations for each value of e. In figure |(||8| we plot q as a function of e for various 
values of the density and different values of k. We see that, as it should be expected for a system undergoing 
a first order phase transition, the curves are smooth for low k and tend to develop a discontinuity for large 
k. For low density the discontinuity occurs at high e, while in is pushed toward small e for high density. 

A different numerical experiment is presented in figure || where we plot the overlap as function of the 
number of Monte Carlo sweeps in a logarithmic scale. Here we let the system evolve at fixed e starting at 
time zero from x = y. Again we observe a behavior compatible with a discontinuity of q as a function of e. 

The same picture can be also observed directly in the behavior of the function gio(r). In figure 10 and 
11 we can see that while gu does not vary too much as a function of e,0 one observes qualitatively different 



regimes in the high and in the low e regime with a discontinuity in <?io(0). 

For high e we observe a strong oscillatory structure in go\, which is very similar to gu, with the only 
difference that the peaks are smoothed and by the presence of the large peak at x = 0. On the contrary 
in the low e regime g\o is very close to one for r > 0.3 and has almost no structure. Notice the striking 
similarity of this picture with the one that we get from the theory. 



4 Notice the small peak at r « v2 for the curves with small e. This could signal some crystallization in the system. In any 
case this effect is present only at high density and is absent in the simulations at lower density. It seems that in the region 
where the overlap becomes quite small, the system crystallize in very long runs, while crystallization is obviously forbidden at 
higher values of the overlap. 
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Figure 9: The overlap as a function of logarithm (in base 2) of the time (i.e. the number of Monte Carlo 
sweeps) starting from y at time and evolving for fixed e. 
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Figure 10: The function <7ii(r) for a system with p = 1.10, N = 1024 and coupling ranging from 0.35 to 
0.6 in steps of 5e = 0.05 The curves do not variate too much with e, except the two with the lowest values 
of e. This could be a spurious effect do to partial crystallization as it is seen from the small peak around 
r = 1.4. 
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Figure 11: The function gio{r) for the same system of figure |H]. We see a net discontinuity in behavior 
from the high to the low e regions. 

6 Conclusions 



In this paper we have shown that constraining the Boltzmann-Gibbs measure to small regions of the 
configuration space we can study the transition from liquid to glass. We have studied in particular the case 
of the simple liquids in the HNC approximation. This model predicts a glassy transition in the lines of the 
Gibbs-Di Marzio scenario. The physical picture got from the present analysis is fully consistent with the one 
obtained from mean- field spin- glass models, although the HNC approach allows to get information on the 
structure function of the supercooled liquid. The HNC approach, which consists in neglecting the so called 
bridge diagrams in the Mayer expansion turns in this way to be equivalent to a kind of mean-field theory, 
where the glass transition is associated to a non-convex effective potential. Many of the prediction got from 
the theory are expected to be valid in real systems. In particular the picture of first order phase transition 
in presence of a coupling, which is associated to the presence of metastable regions in configuration space. 
The picture has been positively tested in numerical simulations on hard sphere systems. 

Growing evidence show that a good starting point for describing the glass transition of supercooled 
liquids is the scenario met in long-range spin glasses with "one step replica symmetry breaking" transition, 
which is the static counterpart of the schematic Mode Coupling Theory. Beyond these models this kind of 
transition pattern is met in mean- field glass models without disorder ||, H [K]] , and simulations of realistic 
glass models confirm many aspect of the picture [11]. Finally the HNC approach of || and with a different 
profile this paper also confirm this scenario. 

We believe that the mean-field theory of the glass transition is now on a firm ground. On one hand we 
have the mode coupling theory that allow to extract the dynamical behavior of supercooled liquids not too 
far from the glass transition. On the other, this theory has been repeatedly observed to be exact in mean- 
field spin glasses, whose study enrich the picture with a static view of the topography of the configuration 
space. In this paper we show that the same picture also hold for the realistic model of liquid obtained by 
the HNC approximation. 

Finally we would like to comment once more on the limitation of the theory. As a genuine mean-field 
theory, the HNC predicts the existence of infinite life metastable states analogously to what happens in long 
range models and in the Mode Coupling approximation. According to this picture the system after having 
frozen into a metastable state at T g , would "follow" it down and back up in temperature in a completely 



reversible manner. While this can be true for short times, irreversible effects are observed on large time 



scales in glasses [21]. The issue is related to the barrier jumping processes that restore ergodicity below T c . 
Unfortunately at present the physics of these dynamical problem is unclear. We hope that insight in this 
problem will come from numerical and theoretical study of disordered models. 
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Appendix 

We show here that the Mean Spherical (MSA) and Percus-Yevick (PY) approximations do not give any glass 
transition. The MSA deals with hard core interactions, for which <f>(x) = oo for \x\ < tq and 4>(x) = 4>i(x) 
for \x\ > tq. In the canonical formalism the exact partition function can be written as a functional integral 
over the density field n(x): 

Z = J Vn(x) exp {/ da; [—/3(fi(x — y)n(x)(n(y) — 6(x — y)) — n(x) logn(x)]} 

S(Jdxn{x)-p) (20) 



In the MSA one substitute the entropic part — J dxn(x) logn(x) in the exponent of ( p0[ ) by the impenetra- 
bility constraint / dxn(x)n(x + y) = if \y\ < r p2| . Expressing the constraint through a set of Lagrange 
multipliers X(x) one gets a Gaussian integral that can be computed exactly. 

In the formalism explained in the text the density variables are replicated, and the constraint takes the 
form 

J dx n a (x)n a (x + y) = for |y| < r$ and all a (21) 

Notice that the replica index in the two densities appearing in (gl]) is the same, as, obviously, there is no 
impenetrability constraint for different replicas like there would be in a multicomponent fluid. Introducing 
now a set Lagrange multipliers A a (x) to enforce (^l|) and exploiting the properties of Gaussian integrals one 
founds that the correlation function 

(n a (x)n b (y)} - (n a (x)) (n b (y)) = A~£{x - y) (22) 

can be found inverting the relation: 

A ab (x) = 5 ab \ a (x) + V ab (x) (23) 

where 

Vab{x) = M x ) 5 ab - {Ml - 5ob)0 b + ($06(1 - 5 0a )p a }ew(x) (24) 

It is clear that for e — > 0, A ab , and therefore its inverse, become diagonal in the replica indices, implying 
9ab{x) = 1 for all x if a ^ b. 

Analogous reasonings apply to the PY approximation (see e.g. ]p!q] ) where one has 

c ab (x) = -g a b(x)[ef 3v ^-l}. (25) 

For e — > c ab becomes diagonal, and the Ornstein-Zernike relation (^) shows that g a b{x) = 1 if a ^ h. 
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